Comparison of mean-field theories for vortices in trapped Bose-Einstein condensates 
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We compute structures of vortex configurations in a har- 
monically trapped Bose-Einstein condensed atom gas within 
three different gapless self-consistent mean-field theories. 
Outside the vortex core region, the density profiles for the 
condensate and the thermal gas are found to differ only by a 
few percent between the Hartree-Fock-Bogoliubov-Popov the- 
ory and two of its recently proposed gapless extensions. In 
the core region, however, the differences in the density pro- 
files are substantial. The structural differences are reflected 
in the energies of the quasiparticle states localized near the 
vortex core. Especially, the predictions for the energy of the 
lowest quasiparticle excitation differ considerably between the 
theoretical models investigated. 
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The landmark experiments to realize Bose-Einstein 
condensation in dilute atomic gases ||l[ have sparked vig- 
orous investigation on the physical properties of these 
novel quantum fluids. Due to the weak interactions, such 
systems are rare examples of interacting quantum flu- 
ids amenable to quantitative microscopic analysis, and 
thus provide unique possibilities to test the fundamental 
principles and theories of many-body quantum physics. 
Theoretical approaches yield several quantities, such as 
density profiles for the condensate and the thermal gas 
component, stability estimates, specific heats, and prop- 
erties of various propagating sound modes, to be com- 
pared with experiments. Experiments also yield detailed 
information on the energies of the individual excitation 
modes of these systems Such information provides 
the most direct and stringent tests for the accuracy of dif- 
ferent theoretical approaches, as compared to the above- 
mentioned "collective" quantities which depend on the 
excitation spectrum as a whole. 

The Bogoliubov equations |^ are a widely used start- 
ing point to compute the excitation spectra for dilute 
Bose-Einstein condensates (BECs). They can be seen 
as eigenmode equations for the condensate described by 
the Gross-Pitaevskii equation Q, neglecting effects of 
the thermal, noncondensed gas component in the sys- 
tem. The Hartree-Fock-Bogoliubov (HFB) theory ||] 
takes self-consistently into account the condensate and 
the thermal gas densities, as well as the lowest order 
anomalous average of the boson field. However, it is 
plagued by an unphysical gap in the excitation spec- 
trum, which violates Goldstone's theorem and invalidates 



its value in predicting the lowest collective mode excita- 
tion frequencies. Goldstone's theorem can be restored 
by neglecting the anomalous average mean field in the 
HFB formalism. This yields the gapless Popov version 
of the HFB theory At low temperatures, predic- 
tions of the Popov approximation (PA) for the lowest 
excitation frequencies of irrotational condensates are in 
good agreement with experimental results, but at tem- 
peratures T ^ Tbec/2 (Tbec denotes the critical temper- 
ature of condensation) the deviations become apparent 
[Q . The main inadequacies of the PA are that it neglects 
the effects of the background gas on atomic collisions 
and the dynamics of the thermal gas component. As an 
improvement to overcome the first limitation within a 
computationally manageable formalism, the so-called Gl 
and G2 approximations have been suggested |^,^ . They 
are gapless mean-field theories which take into account 
effects of the medium on atomic collisions by allowing the 
interaction couplings to depend on the correlation mean 
fields in a self-consistent manner. The two versions are 
based on different approximations for the momentum de- 
pendence of the full many-body T-matrix in the homo- 
geneous gas limit, and their precision for inhomogeneous 
systems remains to be investigated. 

To assess the accuracy of the above-mentioned gapless 
HFB-type approximations, their predictions for the exci- 
tation frequencies of irrotational, harmonically trapped 
atomic BECs have been computed and compared with 
experiments §,|,|l0|-|2l. For temperatures T < Tbec/2, 
the predictions of the PA, Gl and G2 for the lowest ex- 
citation frequencies differ only a few percent HJ^]. For 
higher temperatures, the differences are larger and ex- 
ceed the experimental uncertainty estimates for measure- 
ments, but none of the theories agrees satisfyingly with 
experiments |8|. However, in this temperature range the 
dynamics of the thermal gas component, which these ap- 
proximations do not take into account, is expected to 
have an increasingly important influence on the excita- 
tion eigenmodes. Consequently, results for irrotational 
condensates remain somewhat inconclusive in determin- 
ing the validity of these theories. 

Recently, vortex states in dilute atomic BECs have 
been experimentally realized [|l3|. Furthermore, by ob- 
serving the precession of vortices, the energy of the low- 
est excitation, the so-called lowest core localized state 
(LCLS), has been measured jldj. Interestingly enough, 
the experimental results for this energy agree well with 
the Bogoliubov approximation , while they definitely 
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disagree with the picture given by the self-consistent 
mean-field theories: the latter predict the energy of the 
precession mode to be positive w.r.t. the condensate state 
p6| , pTt , but experiments imply negative energies. We 
suggest that this puzzling fact could be due to incomplete 
thermalization of the (moving) vortex and/or the limita- 
tions of the quasi-stationary, i.e., adiabatic HFB formal- 
ism in describing time-dependent phenomena. The adi- 
abatic approximation essentially fails if the kinetic rates 
of the system exceed frequency separations of the exci- 
tations. Lately, we have shown that the requirement of 
adiabaticity leads to a criterion for the velocity of the 
moving vortex, which is violated in the precession obser- 
vations so far [Q. However, if the precession radii — and 
thus the velocities — of the vortices could be reduced, or 
the physical parameter values appropriately adjusted in 
order for the system to better fulfill the criteria for adia- 
baticity and thermalization, one should be able to mean- 
ingfully compare experimental data with the predictions 
of the self-consistent equilibrium theories for the vortex 
states. 

In this paper, we present results of computations for 
the structures and excitation frequencies of vortex states 
within the Gl and the G2 approximations, and compare 
them with the previously computed predictions of the 
Popov approximation Outside the vortex core 

region, the density profiles for the condensate and the 
thermal gas component are found to differ by only a few 
percent between the PA, Gl and G2. However, in the 
core region the differences are considerably larger. This 
is reflected in substantial differences in the energy of the 
LCLS, which is localized in the core region. 

The gapless HFB-type theories considered in this paper 
can be expressed in the form of the generalized Gross- 
Pitaevskii (GP) equation ||,^ 

[Hoir) + C/c(r)|^(r)p + 2C/e(r)p(r)]^(r) = M</.(r) (1) 

for the condensate wavefunction (t>{r), and the eigenvalue 
equations 

C{r)ug{r) + f/e(r)02(r)z;,(r) = EqUq{r), (2a) 
Cir)vq{r) + i7c(r)0*2(r)u,(r) = -E^v^ir) (2b) 

for the quasiparticle amplitudes Uq{r), Vq{r), and 
eigenenergies Eq. Above, Ho{r) ~ — ?i^V^/2m -I- Vtrap(r) 
is the bare single-particle Hamiltonian for atoms of mass 
m confined by a harmonic trapping potential Vtiap(r) = 
^m{uj^r^ + ijjIz^) expressed in cylindrical coordinates 
r — {r,9,z). Furthermore, /i denotes the chemical po- 
tential, p(r) is the noncondensate density, and £(r) = 
Ho(r) - H + 2C/e(r)|0(r)|2 + 2C/e(r)p(r). The coupling 
functions Uc{r) and Ue{r) specify the approximation: In 
the Popov version Uc{r) = Ue{r) = g, where the coupling 
constant g is related to the s-wave scattering length a by 
g = iTrfi^a/m. In the Gl approximation, Ue{v) = g and 



Uc{r) = g[l + A(r)/(/)^(r)], where A(r) is the anoma- 
lous average of two Bose field operators, and in the G2, 
L/e(r) = Uc{r) = g[l + A(r)/02(r)]. Finally, the selfcon- 
sistency equations 

pir) = + Mr)\')niEq) + \vqir)\% (3) 

A(r) = Y,[2uq{r)v*q{r)n{Eq) + Uq{r)v;{r)], (4) 

where n{Eq) = {e^i^^^'^ — 1)^^ is the Bose distribution 
function (the chemical potential is set equal to the con- 
densate eigenenergy — for the parameter values used in 
the computations, the effect of this approximation is neg- 
ligible), relate the noncondensate density and the anoma- 
lous average to the positive-norm quasiparticle eigenso- 
lutions of equations (H). The expression for the anoma- 
lous average is ultraviolet divergent, and we renormalize 
it by subtracting the last, only implicitly temperature- 
dependent terms from the sum of equation (^. This 
renormalization scheme neglects the difference between 
the two- and the many-body T-matrices, but for dilute 
gases the corrections are vanishingly small . 

In this paper we consider axisymmetric single-quantum 
vortex states of the form 0(r) = (/)(r)e*^, thus neglecting 
possible bending of the vortex [ pO| . Furthermore, we re- 
strict to the case lOz — 0, which implies the system to be 
homogeneous in the axial direction. However, the qual- 
itative results concerning especially the substantial rela- 
tive differences between the approximations in the lowest 
excitation energy are expected to remain valid also for 
systems with > 0. 

In accordance with the translational symmetry in the 
axial direction, we impose periodic boundary conditions 
at z = ±L/2. By substituting the Ansatz 



uq{v) = M^(r)e"'-(2^/^)^+»(9«+i)«, 
Vq{v) = t;^(r-)e^9^(2Vi)^+%e-i)e^ 



(5a) 
(5b) 



where qe and are integer quantum numbers, equa- 
tions (H) reduce to radial equations. We discretize 
them with a high-order finite-difference method and 
solve the consequent matrix eigenvalue problem by us- 
ing an Arnoldi method implemented in the ARPACK 
( http : / / www, caam . rice . cdu/software /ARPACK / ) sub- 
routine libraries [[17[. The nonlinear Gross-Pitaevskii 
equation is solved using a real-space discretization and 
iterative relaxation methods. Finally, self-consistency is 
achieved by iteration: Using the solutions of equations 
(0) and (^) corresponding to a given total number of par- 
ticles, new mean-field potentials are computed from the 
self-consistency equations, and the procedure is repeated 
until convergence. We solve quasiparticle states explic- 
itly up to an energy i?c, and compute the contribution to 
the mean-field potentials from the states above E^ by us- 
ing a semiclassical local density approximation described 
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in reference |ll|]. This method allows one to use rather 
small values for Ec with excellent accuracy, thus essen- 
tially improving the computational efficiency. In order to 
stabilize the iteration, we use underrelaxation in updat- 
ing the mean-field potentials. 

To facilitate comparison with previously presented re- 
sults for the Popov approximation, the physical parame- 
ter values for the gas and the trap were chosen to be the 
same as those in reference jl^]. We modelled a sodium 
gas with the atomic mass m = 3.81 x 10~^^ kg and the 
scattering length a = 2.75 nm in a trap with the radial 
frequency i^r = Wr/27r = 200 Hz. The density of the gas 
was determined by treating N = 2x 10^ atoms per length 
L = 10 /im in the axial direction. Altogether, these val- 
ues yield the condensation temperature Tbec ~ 0.8 /iK. 

Figures and ^ present results of our computations for 
axisymmetric single-quantum vortex states. The density 
profiles for the condensate, the noncondensate and the 
anomalous average are displayed in figure |l| at temper- 
atures r = 50 nK and 400 nK. Outside the vortex core 
region, the differences in the density profiles between the 
PA, Gl and G2 are at most a few percent at temperatures 
T < Tbec/2. In the core region, however, the differences 
are considerably larger. The total density of the gas is 
approximately 20% larger on the vortex axis in the G2 
than within the PA. This squeezing behavior is associated 
to the "softening" of the repulsive effective interaction in 
the core region due to many-body effects ^ . 

The differences in the core densities between the var- 
ious approximations also suggest differences in the ener- 
gies of the quasiparticle excitations localized in the core 
region. Figure ^ displays the energies of three such states, 
the lowest excitations with angular momentum quantum 
numbers qe = —1,0, and 1, as functions of tempera- 
ture. For the Gl and the G2, the increased core densities 
are compensated by smaller effective couplings, and the 
shifts in the excitation energies are generally only a few 
percent for temperatures T ^ Tbec/2; at higher tempera- 
tures the softening effect of the interaction becomes more 
pronounced |^,^, also increasing the shifts in the exci- 
tations. However, the lowest Kelvin mode (consisting of 
the lowest qg — —1 excitations) state, the LCLS, is espe- 
cially sensitive to the structure of the core region. The 
differences in the energies of the LCLS between the ap- 
proximations are 25-40% even at temperatures for which 
the predictions of the PA for the excitation frequencies 
of irrotational condensates differ by less than 5% from 
the experimental data. In addition, the temperature de- 
pendence of the LCLS is found to be much stronger than 
for the other states. In fact, the lowest excitation en- 
ergy vanishes in the zero-temperature limit for all the 
approximations [ p7| ; the remainder of the spectrum is 
essentially temperature- independent, except in the vicin- 
ity of Tbec- The State displayed in figure |^ with qe = I 
is the Kohn mode, which should have the exact energy 
E = hujr according to Kohn's theorem for parabolic con- 



finement |gl|. Kohn's theorem is satisfied to an accuracy 
of 1-6% for all the approximations, suggesting that dy- 
namical effects of the thermal gas component are small 
in the temperature range studied. 

In conclusion, we argue that future measurements of 
the lowest excitation frequencies of the vortex states 
could provide stringent tests for the validity of the mean- 
field theories considered. Especially, they could be used 
to estimate the degree to which the approximations for 
the many-body T-matrices based on the homogeneous 
limit remain valid for highly inhomogeneous systems. 
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FIG. 1. Density profiles of the vortex state for the conden- 
sate (|<^|^), thermal gas component (p) and anomalous average 
(|A|) in the PA (solid), Gl (dashed) and G2 (dashed-dotted) 
at temperatures (a) T = 400 nK and (b) T = 50 nK. Axes 
for the values of p and |A| are on the left-hand sides, and 
for \(t>'\^ on the right-hand sides. Figures (c) and (d) display 
the vortex core region, where the differences in the density 
profiles between the approximations are substantially larger 
than farther from the vortex axis. 



FIG. 2. (a) Energies of the lowest excitation modes with 
angular momentum quantum numbers qe — —1,0, and 1 in 
the PA (solid), Gl (dashed) and G2 (dashed-dotted) as func- 
tions of temperature. The qe = ~1 state is the so-called 
lowest core localized state (LCLS), and the qg = 1 excitation 
is the Kohn mode, (b) Temperature dependence of the energy 
of the LCLS within the PA, Gl and G2. Note the substantial 
relative differences between the theories in this energy. 
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